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Abstract 

As global warming accelerates the melting of Arctic sea ice, polar bears (Ursus 
maritimus) must adapt to a rapidly changing landscape. This process will neces- 
sarily alter the species distribution together with population dynamics and 
structure. Detailed knowledge of these changes is crucial to delineating conser- 
vation priorities. Here, we sampled 361 polar bears from across the center of 
the Canadian Arctic Archipelago spanning the Gulf of Boothia (GB) and 
M'Clintock Channel (MC). We use DNA microsatellites and mitochondrial 
control region sequences to quantify genetic differentiation, estimate gene flow, 
and infer population history. Two populations, roughly coincident with GB and 
MC, are significantly differentiated at both nuclear (F sr = 0.01) and mitochon- 
drial ($ ST = 0.47; F sr = 0.29) loci, allowing Bayesian clustering analyses to 
assign individuals to either group. Our data imply that the causes of the mito- 
chondrial and nuclear genetic patterns differ. Analysis of mtDNA reveals the 
matrilineal structure dates at least to the Holocene, and is common to individ- 
uals throughout the species' range. These mtDNA differences probably reflect 
both genetic drift and historical colonization dynamics. In contrast, the differ- 
entiation inferred from microsatellites is only on the scale of hundreds of years, 
possibly reflecting contemporary impediments to gene flow. Taken together, 
our data suggest that gene flow is insufficient to homogenize the GB and MC 
populations and support the designation of GB and MC as separate polar bear 
conservation units. Our study also provide a striking example of how nuclear 
DNA and mtDNA capture different aspects of a species demographic history. 



doi: 10.1002/ece3.662 

Introduction 

Climate change is expected to have a profound impact on 
polar bears (Ursus maritimus), affecting the spatial and 



seasonal distribution of the sea ice that they rely on to 
travel and hunt (Stirling and Derocher 1993; Derocher 
et al. 2004). During the 21st century the availability of 
suitable polar bear habitat is predicted to decline, with a 
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potential negative effect on the species abundance (Amst- 
rup et al. 2008; Durner et al. 2009). Adaptation to this 
rapidly changing environment, should it occur, will neces- 
sarily involve changes to the spatial distribution of 
individuals. This may result in range shifting or local 
extinction of some populations. The design of cogent and 
effective conservation strategies requires a detailed under- 
standing of how the population dynamics are changing. 
This is of even greater importance if traditional harvesting 
practices of the species are to continue on a presumed 
sustainable basis (Lee and Taylor 1994; Schliebe et al. 
2008). 

Currently, 19 distinct polar bear populations are recog- 
nized worldwide (Obbard et al. 2010; see Paetkau et al. 
1999 for a detailed map of the geographic distributions of 
these populations), established on the basis of land barri- 
ers and movement patterns inferred through mark-recap- 
ture and radiotelemetry data (Taylor and Lee 1995; 
Bethke et al. 1996; Taylor et al. 2001). These 19 popula- 
tions constitute the management units that are used to 
delineate global conservation strategies and establish har- 
vest quotas. Evidence from genetic data also supports the 
designation of the majority of these populations (Paetkau 
et al. 1995, 1999). In the most geographically comprehen- 
sive study on polar bear genetics performed to date, Pae- 
tkau et al. (1999) used 16 microsatellite loci to survey 16 
of the 19 recognized polar bear populations. Although 
there were significant genetic differences for most pairwise 
population comparisons, these were generally small with 
no striking discontinuities across the range. Populations 



clustered into four groups suggesting some overall genetic 
structure. Paetkau et al. (1999) interpreted these differ- 
ences as the product of contemporary movement patterns, 
possibly reflecting the configuration of land masses in 
conjunction with the seasonal distribution of sea ice and 
how this enables access to seals. 

The study by Paetkau et al. (1999) focused on obtain- 
ing samples from nearly all polar bear populations, thus 
providing insight into the genetic structure across the 
species' range. An alternative strategy, which targets a 
small number of populations with continuous and inten- 
sive sampling (e.g., Cronin et al. 2006; Crompton et al. 
2008; Zeyl et al. 2009), allows us to characterize genetic 
structure in greater detail and assess how these groups 
are or have been connected and perhaps how they arose. 
Here, we focus on two previously designated polar bear 
populations geographically located in the center of the 
Canadian Arctic Archipelago: the Gulf of Boothia (GB), 
which encompasses an area of -170,000 km 2 ; and the 
adjacent M'Clintock Channel (MC), that extends over 
~ 500,000 km 2 . These areas are separated by a land bar- 
rier (the Boothia Peninsula) for the majority of their 
length and connected by water only through the narrow 
Bellot Strait, an ~25 km long passage with fast tidal cur- 
rents that prevent the formation of stable pack ice 
(Fig. 1). While the GB polar bear population has been 
suggested to be demographically stable, the MC popula- 
tion is increasing in size due to management measures 
that have allowed recovery from intense hunting pressure 
(Obbard et al. 2010). These populations are situated 



Figure 1. Map of the study area showing the 
locations where each of the 718 polar bears 
where captured. The dotted line represents the 
border between the Gulf of Boothia and 
M'Clintock Channel populations, designated 
through previous studies of polar bear 
movement patterns. 




200 km 
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within the Northwest Passage, a new navigation route 
that will increasingly facilitate shipping across the Arctic 
seas along the top of North America as sea ice continues 
to melt (Kerr 2002). This expanding development of the 
Arctic puts the GB and MC polar bear populations at 
risk of being significantly disturbed by anthropogenic 
stressors (e.g., increased pollution or direct interactions 
with humans) in the near future (Amstrup et al. 2008; 
Obbard et al. 2010). In our study, we use mitochondrial 
and nuclear DNA sequences combined with extensive sam- 
pling to address the following questions: 

(1) To what extent are the GB and MC populations 
genetically differentiated? If significant genetic differ- 
ences do exist, over what time frame did they arise 
and what might have been the cause? 

(2) Does gene flow occur between regions and if so, is 
there evidence of sex-biased gene flow? Do diagnosed 
patterns of gene flow match our knowledge of move- 
ment patterns deduced from mark-recapture data? 

Material and Methods 
Sampling and data set 

Samples for genetic analyses were collected by the Nunavut 
Department of Environment during surveys carried out in 
the spring (March-June) between 1998 and 2000 in two 
previously defined polar bear populations: GB and MC 
(Fig. 1). Part of this region is characterized by perennial 
ice where polar bears remain throughout the year (Amst- 
rup et al. 2008; see Barber and Iacozza (2004) for a 
detailed analysis of the annual sea ice patterns in the study 
area). After bears were immobilized, capture site and sex 
were recorded, age was estimated and a tissue sample (a 
small disk of skin obtained from ear tagging) was taken 
for DNA analysis. For details on sampling procedures, see 
Taylor et al. (2006, 2009). In total our data set contains 
718 samples, including mothers with their cub(s) and soli- 
tary males; some individuals were captured up to three 
times (during each of the 3 years of fieldwork). For genetic 
analyses we used 361 unique male and female adult and 
subadult individuals (GB, N = 289; MC, N = 72), discard- 
ing genetic data from cubs, yearlings and second-year 
bears that were related to the sampled mother. 

We used mark- recapture data to compare movement 
patterns between sexes and populations. Movement was 
estimated by calculating the straight-line distance 
between sites where each individual was captured in dif- 
ferent years using the great circle distance calculator in 
Google Maps. When a bear was captured in three differ- 
ent years, we only used the maximum distance between 
capture sites to avoid pseudo-replication. Three females 
moved from MC to GB and could not be a priori 
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assigned to either population (and thus be confidently 
used to compare movement patterns between GB and 
MC). We analyzed the data in two ways, first excluding 
these three females and then again considering them part 
of the population where they were initially captured 
(MC). In total we obtained 56 independent observations 
of distances traveled by individual bears during the study 
period. We compared movement patterns between sexes 
and populations by performing a two factor analysis of 
variance (ANOVA) in JMP version 10 (SAS Institute 
Inc., Cary, NC). 

Microsatellite amplification and genotyping 

Genomic DNA was extracted from a small disk of skin 
using the QIAGEN DNA extraction kit (QIAGEN, Missis- 
sauga, Canada) following the manufacturer's procedures. 
We genotyped 531 individuals (including 361 adults/subad- 
ults, and an additional 170 cubs, yearlings and second-year 
bears for maternity analysis) for nine previously published 
dinucleotide microsatellite loci (G1A, G1D, G10B, and G10L: 
Paetkau and Strobeck 1994; G10M, G10P, and G10X: 
Paetkau et al. 1995; MU59: Taberlet et al. 1997; G10H: 
Paetkau et al. 1998), following protocols outlined in Saun- 
ders (2005). Random subsets of individuals were genotyped 
more than once to confirm repeatability of results. We used 
the family groups in our data set to evaluate possible geno- 
typing errors by performing a maternity analysis with 
Cervus version 3.0.3 (Kalinowski et al. 2007). This data set 
includes 170 cubs belonging to 108 known mothers, thus 
totaling 2502 genotypes across our nine loci. The maternity 
analysis revealed 24 single locus genotypes (0.96% averaged 
across all genotypes or 2.14% averaged across loci) that did 
not match between cubs and their mothers. One mismatch 
could have been caused by null alleles, but the remaining 
23 could be due to either mutation or genotyping error. 
Regardless, we were confident that this low rate would have 
negligible impact on our results. Adult and subadult indi- 
viduals were used to test for Hardy-Weinberg equilibrium 
and linkage disequilibrium in Arlequin version 3.5.1.3 
(Excoffier and Lischer 2010); we excluded all cubs to reduce 
biases due to inclusion of related individuals. The nine loci 
show no deviations from Hardy-Weinberg expectations, 
both when individuals from GB and MC were analyzed 
together and separately. When GB and MC individuals 
were analyzed together, the locus pairs G10M/MU59 and 
G10B/G10H showed significant linkage desequilibrium 
after sequential Bonferroni correction (Rice 1989). When 
GB and MC samples were analyzed separately, only G10B 
and G10H showed significant linkage disequilibrium in the 
GB population, perhaps a consequence of admixture 
between genetically distinct individuals occurring within 
GB. 
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Mitochondrial control region amplification 
and sequencing 

We selected a subset of 86 subadult and adult individuals 
evenly distributed throughout the study area (GB, 
N = 48; MC, N = 38) for the amplification and sequenc- 
ing of an -472 base pair (bp) fragment of the mitochon- 
drial control region (CR; positions 16,545-17,018 in the 
polar bear mitochondrial genome reference sequence 
AF303111 deposited in GenBank). Amplification was car- 
ried out using primers PBCR1 (5'-AGCTCCACTACCAG 
CACCC-3') and PBCR 4 ( 5' - AAATGC ATGAC ACCACA 
GTTATGTGTGATC-3') (Weber et al. unpubl. ms.). Poly- 
merase chain reaction (PCR) cocktails contained 2 /(L of 
genomic DNA, 2.5 mmol/L MgCl 2 , 0.1 mmol/L dNTPs, 
0.1 /.imol/L of each primer, and 0.75 U taq DNA poly- 
merase (Vivantis Technologies, Malaysia) in a final vol- 
ume of 25 /jL PCR ViBuffer A (Vivantis; 50 mmol/L KC1, 
10 mmol/L Tris-HCl pH 9.1 and 0.01% Triton X-100). 
The thermocycling profile included an initial denaturation 
step at 94°C for 5 min; followed by 34 cycles of 94°C for 
30 sec, 50°C for 30 sec, and 72°C for 45 sec; and a final 
extension at 72°C for 5 min. An aliquot of the PCR prod- 
ucts was visualized on an ethidium bromide-stained 2% 
agarose gel and successful reactions were purified using 
the QIAquick PCR purification Kit (QIAGEN). Sequenc- 
ing was carried out in both directions using the above 
mentioned primers at the London Regional Genomics 
Centre (London, Ontario, Canada) using an Applied 
Biosystems 3730 Analyzer. All sequences were deposited 
in GenBank (Accession Numbers KF192517 - KF192602). 

Measures of genetic distance 

Genetic differentiation between populations was estimated 
with Arlequin version 3.5.1.3, using DNA microsatellite 
and CR data to calculate different F-statistics or their ana- 
logs. For DNA microsatellites, we calculated F ST values; for 
the CR data we estimated differentiation using haplotype 
frequencies alone (jF st ), or taking into account the differ- 
ence among haplotypes together with their frequencies 
($ ST ). Statistical significance was assessed through 1000 
random permutations. We calculated the genetic diversity 
within populations by estimating p-distances in Arlequin 
and _F IS values in FSTAT version 2.9.3.2 (Goudet 1995), for 
sequence and microsatellite data, respectively. We used the 
program Populations version 1.2.32 (Langella 1999) and 
microsatellite data to calculate D sw pairwise genetic 
distances (Shriver et al. 1995; a measure of genetic distance 
developed for microsatellite loci) between all individuals. 
The program GenAlEx version 6 (Peakall and Smouse 
2006) was used to calculate geographic distances between 
individuals from geographic coordinates of capture sites. 



We tested for the effect of individual-based isolation by 
distance (Rousset 2000) by performing a Mantel test (Man- 
tel 1967) in GenAlEx, with the null hypothesis of indepen- 
dence between the genetic and geographic distance 
matrices, and assessing significance through 999 random 
permutations. We analyzed the data separately for males, 
females, and both sexes pooled, and performed these analy- 
ses for GB and MC separately, and for both regions pooled. 
Because nine comparisons were carried out sequential Bon- 
ferroni corrections were applied to minimize type I errors. 

Bayesian clustering analyses 

We assessed population structure using individual geno- 
types for the nine microsatellite loci and two programs 
that implement Bayesian clustering algorithms: Structure 
version 2.3.4 (Pritchard et al. 2000) and Geneland version 
4.0.3 (Guillot et al. 2005). Analyses were conducted for 
males and females separately and pooled. Structure was 
run using the admixture ancestry model, correlated allele 
frequencies and both with and without LOCPRIOR (i.e., a 
prior that indicates the population sampling origin). We 
explored values of K = 1 through 4 (two more than the a 
priori reported number of populations) with ten itera- 
tions per value of K each with 2,500,000 generations, 
discarding the initial 500,000 as burn-in. The most likely 
K value was determined using the Evanno et al. (2005) 
method implemented in Structure Harvester version 0.6.93 
(Earl and vonHoldt 2012), and by inspecting individual 
assignment patterns (specifically when comparing between 
results from K = 1 and K = 2, for which the Evanno 
et al. method cannot be applied). 

We conducted two independent Geneland runs per 
group analyzed (males, females, and sexes pooled), each 
consisting of 2,500,000 iterations and a thinning of 100. 
The number of genetic populations was set to K = 2 
based on results of Structure runs, using correlated allele 
frequencies and the spatial model. After discarding 5000 
steps (25%) as burn-in, we divided the study area into 
100 quadrats (20 by 5, -1° longitude by 1° latitude; 
Fig. 2D-E), to obtain a map of population membership 
and F ST values between the two genetic clusters. Similar 
results were reached for replicate runs. 

Network and phylogenetic analyses 

We aligned mitochondrial sequences using BioEdit version 
7.0.9.0 (Hall 1999) and explored the relationship between 
CR haplotypes using a statistical parsimony network 
analysis in TCS version 1.21 (Clement et al. 2000). To 
compare the variation found in CR sequences obtained 
from individuals belonging to GB and MC to the overall 
phylogeographic structure observed in the species, we 
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Figure 2. Results from Bayesian clustering analyses performed in Structure and Geneland using data from nine DNA microsatellite loci. Structure 
bar plots representing assignment of adult and subadult individuals (A), females (B), or males (C) from GB and MC to one of two genetic 
populations (represented by light and dark gray bars). Admixed individuals are represented by bars with varying proportions of light and dark gray 
(scale on the y-axis). Maps of the study area divided in quadrats of -1° longitude by 1° latitude that were assigned to one of two genetic 
populations (light or dark gray), according to the genotypes of sampled individuals (D: adult and subadults; E: females; F: males). Black squares 
represent sampled individuals; thus assignment of quadrats where individuals have not been sampled represent extrapolations and must be 
interpreted with caution. The dotted line shows the border between the GB and MC populations. F ST values represent comparisons between 
individuals from GB and MC (denoted by F ST [pop]) calculated in Arlequin or degree of divergence between genetic clusters (denoted by F S t 
[cluster]) calculated in Geneland. 



downloaded published CR sequences (Edwards et al. 
2011) from Genbank. In total, we incorporated 45 
sequences obtained from modern individuals sampled 
from across the Arctic (from 1883 to present, accession 
numbers: GU573485; GU573490; and sequences with 
locality information with numbers between JF900105- 
153). Sequences were aligned to our 86 CR sequences 
from GB and MC, creating a combined data set of 131 
sequences from individuals across the species' range with 
which we constructed a new statistical parsimony network. 
We also used this combined CR data set to build a Bayes- 
ian phylogenetic tree using MrBayes version 3.2.1 (Huel- 
senbeck and Ronquist 2001; Ronquist and Huelsenbeck 
2003). The tree was rooted using a brown bear CR 
sequence {Ursus arctos: GenBank accession number 
EF033706). The HKY (Hasegawa et al. 1985) +G model of 
nucleotide evolution best fit the data according to an anal- 
ysis performed in jModeltest version 0.1.1 (Posada 2008). 
Two simultaneous Bayesian searches using four incremen- 
tally heated Markov chains and default priors for all 
parameters were run for six million generations. At this 
point the standard deviation of split frequencies was <0.01, 
suggesting the analysis had converged. Trees were sampled 
every 100 generations, and we discarded the first 25% as 
burn-in. All parameters had a potential scale reduction fac- 
tor (Gelman and Rubin 1992) that was close to one, indi- 



cating that we had adequately sampled the posterior 
distribution. We also assessed convergence using the 
"cumulative" and "compare" functions implemented in 
the software AWTY (Wilgenbusch et al. 2004), confirming 
that runs had reached stationarity. Finally, we obtained a 
50% majority rule consensus from the retained, combined 
posterior tree distribution of both runs. 

Bayesian coalescent simulations 

We used the program IMa2 (Hey 2010), that implements 
the isolation with migration (IM) model, to estimate gene 
flow, effective population sizes, and divergence time 
between GB and MC. We ran the program using either 
mitochondrial DNA (mtDNA, CR sequences) or nuclear 
DNA (nuDNA, nine DNA microsatellite loci) separately, 
applying the HKY model to the former and the stepwise 
mutation model to the latter (Ohta and Kimura 1973). In 
both cases the IM model was simplified to include only one 
(bidirectional) migration rate parameter, thus reducing 
the total number of parameters estimated by the model. 
Positions in the CR sequence alignment with ambiguous 
bases or missing data were ignored by IMa2; thus the 
analysis was conducted using 271 bp of the 472 bp ampli- 
fied originally. Runs in M mode showed adequate mixing 
with 100-120 chains and a burn-in of 250,000-500,000 
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generations for the CR and microsatellite data, respectively. 
We implemented the geometric heating model and ran the 
program at least four times, with different random seeds, 
until a minimum of 100,000 genealogies were saved. Joint- 
posterior densities of model parameters were finally 
estimated in L mode. Estimations of population migration 
rates (2Nm), the effective number of migrants per genera- 
tion, were obtained from the migration (m//z) and 9 (4N/J.) 
parameters calculated with IMa2. By calculating 2Nm = 
AN [i x 1/2 x m//.i, we estimated 2Nm independent of the 
mutation rate (Hey and Nielsen 2004). Because the model 
was simplified to include only one bidirectional migration 
parameter, we used the average 9 between GB and MC to 
calculate 2Nm. Finally, effective population sizes for GB 
and MC were also estimated using DNA microsatellite data 
and an approximate Bayesian computation framework in 
ONeSAMP (Tallmon et al. 2008), exploring parameter val- 
ues up to each population's census size. 

Results 

Movement patterns inferred through 
mark-recapture data 

Based on our mark-recapture data, the average distance 
traveled was slightly larger for males than for females 
(71 ± 60 km vs. 57 ± 61 km), and for MC than for GB 
individuals (89 ± 90 km vs. 54 ± 45 km), although these 
differences were not statistically significant. The maximum 
distance traveled was larger for females (317 km) than for 
males (194 km), and for MC (317 km) than from GB 
(185 km). In three instances, females that were initially 
observed within the MC boundaries were subsequently cap- 
tured in the area encompassed by GB (having traveled 86, 
137, or 230 km). When these three females were included in 
the population from which they had initially been captured 
(MC), MC individuals were found to have traveled larger 
distances compared to those of GB (101 ± 88 km vs. 
54 ± 45 km), and this difference was statistically significant 
(two factor ANOVA: F h 52 = 6.310, P = 0.015). 

Nuclear differentiation between individuals 
from GB and MC 

Bayesian clustering analyses performed using data from 
nine microsatellite loci suggest GB and MC polar bear 
populations are genetically differentiated (Fig. 2). The 
most likely number of genetic populations (K) identified 
by Structure was two, corresponding to differentiation 
between the GB and MC populations, although 
individuals show varying levels of admixture (Fig. 2A). 
Using the Evanno et al. (2005) method, we did not find 
support for values of K greater than 2, and likelihood 



scores were slightly larger for K = 1 than K = 2 (see 
Appendix SI). However, for K=2 Structure assigned 
-90% of the membership of GB and -80% of MC into 
two distinct clusters (Fig. 2 A and Appendix SI), suggest- 
ing genetic structure can be detected between these popu- 
lations. The differentiation was weak, although statistically 
significant (F S t = 0.01, P < 0.0001), causing Structure to 
fail in assigning individuals to either genetic population 
without incorporating sampling locality as a prior in the 
Bayesian analysis. Using information from individual 
genotypes and sampling coordinates, Geneland assigned 
most quadrats in the GB to the GB genetic cluster, and 
most MC quadrats to the MC cluster, with an F S t value 
of 0.013 between clusters (Fig. 2D). However, as a conse- 
quence of the apparent admixed ancestry of some individ- 
uals, a few quadrats within the GB were more closely 
allied to the MC cluster, and vice versa (Fig. 2D). We 
found no isolation by distance, either when all samples 
were analyzed together or when GB and MC samples 
were analyzed separately. Note that not every quadrat was 
sampled (Fig. 2D-F), and thus the assignment of such 
"empty" quadrats represents an extrapolation that must 
be interpreted with caution. Moreover, the genetic iden- 
tity of quadrats within populations must be interpreted as 
transitory. Geneland performs these assignments using the 
coordinates where the individual was sampled and it is 
likely that polar bears move freely within GB and MC. 
Moreover, this analysis was conducted pooling individuals 
from three different years and our mark-recapture data 
suggest that average movements on the order of 
50-70 km between years are common. Thus our Geneland 
results are approximate and it is likely that the genetic 
identity of quadrats within GB and MC changes as bears 
move on the landscape. Finally, the Fi S value was slightly 
larger in MC than in GB (0.012 vs. -0.001). Altogether, 
these results suggest nuclear genetic differentiation 
between GB and MC with a low level of admixture. 

We tested for sex-specific patterns in DNA microsatel- 
lites by performing the same analyses as above for each 
sex separately. Again Structure assigned individuals to two 
genetic clusters corresponding to GB and MC (when 
capture location was used as a prior); however, differenti- 
ation was higher in females than males (Fig. 2B vs. 2C). 
-Fst values between GB and MC individuals and between 
diagnosed genetic clusters (i.e., the differentiation between 
light and dark gray quadrats identified by Geneland) were 
both marginally higher for females than for males (0.01 
vs. 0.009 and 0.016 vs. 0.01, respectively). As before, for 
both sexes Geneland assigned most quadrats in GB to a 
different genetic cluster than those in MC and vice versa 
(Fig. 2E and F). Again, we found no significant isolation 
by distance for males or females overall or analyzing GB 
and MC individuals separately. Although sample sizes are 



2013 The Authors. Ecology and Evolution published by John Wiley & Sons Ltd. 



3157 



Polar Bear Genetics in the Canadian Arctic 



L. Campagna ef al. 



lower for males compromising statistical power, these 
results are consistent with males exhibiting higher 
dispersal between populations. 

Mitochondrial differentiation between 
individuals from GB and MC 

The CR statistical parsimony network analysis implies 
genealogical distinction between GB and MC (Fig. 3). The 
network analysis revealed 10 CR haplotypes (labeled a-j) 
among 86 individuals that differed from each other by up 
to nine mutational steps (Fig. 3A). Haplotypes and their 
frequencies differed significantly between GB and MC 
($ ST = 0.47 and F ST = 0.29, P < 0.001), with haplotype e 
being most common among GB individuals and haplotype 
a most abundant among MC individuals (Fig. 3A). The 
genetic diversity within both populations was similar (aver- 
age p-distance values of 2.1 ± 2.6 in MC vs. 1.2 ± 1.9 in 
GB). Figure 3B maps the location where polar bears 
carrying the two most frequent haplotypes (a and e) 
were captured; again showing how individuals carrying the 
e haplotype are more common in GB and those with the 
a haplotype are mainly found in MC. 

The divergence between individuals from GB and MC is 
comparable to the divergence among CR sequences from 
individuals sampled across the species' range (Fig. 4). 
When the CR data were reanalyzed together with 45 new 
sequences obtained from sites shown in Figure 4C, the 
statistical parsimony analysis (Fig. 4A) identified 22 CR 
haplotypes separated by up to 11 mutational steps. Thus 
the diversity in CR sequences observed at the scale of 
the central portion of the Canadian Arctic Archipelago is 
similar to the species CR genetic variation as a whole 
(compare the haplotype networks in Figs. 3A, 4A). Neither 
the network analysis nor the Bayesian tree (Fig. 4B) reveals 
a completely resolved genealogical structure. Haplotypes 
from the same or adjacent localities are often scattered 
throughout the parsimony network in Figure 4A. Although 
support for nodes in the Bayesian tree was low (generally 
<0.95) suggesting more data are needed to obtain a fully 
resolved topology, individuals from the same locality can 
be found in nearly every clade (Fig. 4B). The only highly 
supported clade within the species (posterior probability 
of 1.00) includes individuals from six of the nine localities 
sampled (Fig. 4B). 

Bayesian coalescent simulations 

Estimates of the demographic parameters from the IM 
model differed between our mtDNA and nuDNA analyses 
(Table 1). Analyses of the mitochondrial data suggested 
that the effective population size is approximately fivefold 
larger in MC than in GB, while the DNA microsatellite data 



implied a threefold size difference in the opposite direction 
(Table 1). The latter trend was also obtained when effective 
population size was estimated using DNA microsatellite 
data and an Approximate Bayesian computation frame- 
work in ONeSAMP (average and 95% confidence interval 
limits: GB, 421 [371-555] individuals vs. MC, 80 [66-111] 
individuals). IMa2 analysis of mtDNA suggests that there 
have been significant levels of gene flow since the establish- 
ment of the GB and MC populations (95% highest poster- 
ior density interval - 95% HPD - does not overlap with 
zero), whereas the converse is true when the analysis was 
performed using nuclear data (Table 1). Although esti- 
mates of splitting time from the IM model must be inter- 
preted with caution due to the lack of species-specific 
molecular clock calibrations, we can use calibrations from 
other large mammal species. Using radiocarbon-dated sam- 
ples and DNA sequences, Saarma et al. (2007) estimated 
the brown bear ( 17. arctos) CR to diverge at a rate of 29.8% 
(confidence interval: 13.3-47.6%) per million years. Huang 
et al. (2002) utilized known human pedigrees to infer a 
mutation rate of 1 x 10~ 4 per generation for dinucleotide 
microsatellites in Homo sapiens. Assuming that these muta- 
tion rates are comparable to those in U. maritimus, and 
using a generation time of 5 years (approximated through 
measurements of age of first reproduction of females in the 
study area, Taylor et al. 2006, 2009), we obtain an approxi- 
mately 40-fold difference in the estimation of the splitting 
time model parameter between populations for the two 
classes of markers. Our analyses suggest differences in 
mtDNA haplotypes began accumulating in the Holocene 
(ca. 3500 years before present), while nuDNA implies GB 
and MC began to diverge very recently, ca. 85 years ago 
(95% HPD: 14-350 years ago). Posterior density curves for 
parameters estimated in the IMa2 analyses and details for 
the calculations regarding the splitting time parameter are 
presented in Appendix S2, S3, respectively. An unlikely 
high divergence rate of ca. 1200% per million years or an 
extremely low mutation rate of 2.5 x 10~ 6 per generation 
(for the CR and microsatellites, respectively) would be 
required for time estimations to be congruent (Appendix 
S3). The posterior density curve for the r/i parameter esti- 
mated using mtDNA shows a sharp peak but does not 
return to zero at the upper bound of the prior (Appendix 
S2). We interpret this as representing a range of times for 
the establishment of the matrilineal structure that are 
equally probable and thus our estimate (obtained using the 
peak value) is conservatively recent. Estimates of effective 
population sizes and gene flow were similar when male and 
female microsatellite data were analyzed separately; how- 
ever, values of tfi differed (Table 1). While the 95% HPD 
interval for t\i in males overlaps with zero, the estimate for 
females is similar to the value obtained when the data were 
pooled (Table 1). 
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Figure 3. Statistical parsimony network 
analysis performed with CR sequences from 86 
polar bear individuals. (A) Maximum parsimony 
network showing 95% probability linkages 
among 10 CR haplotypes. Lines represent 
mutational changes with black squares 
indicating hypothetical haplotypes. The area of 
each circle is proportional to the number of 
individuals with that haplotype: the smallest 
circles are haplotypes found in one individual 
while the largest represents 33 individuals. 
Haplotypes are coded in dark or light gray 
following the populations where they were 
sampled (GB or MC, respectively). (B) Map 
representing the sampling location of 
individuals carrying the two most common 
haplotypes (a and e) which differ in frequency 
across populations. 
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Discussion 

Our analyses of mitochondrial and nuclear sequence data 
suggest that polar bears from the GB and MC do indeed 
comprise genetically differentiated populations consistent 
with previous designations (Taylor and Lee 1995; Bethke 
et al. 1996; Paetkau et al. 1999; Taylor et al. 2001), albeit 
at low levels that imply ongoing gene flow and/or recent 
common ancestry. Although analyses of both mtDNA and 
nuDNA support this assertion, these distinct classes of 
molecular markers reveal signals of divergence at different 
temporal and spatial scales. We found significant differ- 
ences in CR haplotype frequencies between the GB 
and MC populations, and mtDNA structure common to 
individuals from across the species' range (Fig. 4) that 
most likely originated elsewhere before or during the 
Holocene. In contrast, the differentiation in DNA micro- 
satellites dates to modern times, possibly reflecting con- 



temporary movement patterns, changes in the species 
distribution, and impediments to gene flow that have 
resulted from recent changes to the Arctic landscape (e.g., 
changes in the ice conditions that influence hunting 
opportunities). Below we discuss our findings in the con- 
text of the questions that motivated this study. 

Population differentiation across multiple 
temporal and spatial scales 

Paetkau et al. (1999) found the GB and MC populations 
to cluster with others from the Canadian Arctic Archipel- 
ago (Lancaster Sound and Viscount Melville Sound) and 
the two adjacent populations to the east (Baffin Bay and 
Kane Basin). Paetkau et al. (1999) reported an _F S t value 
of 0.011 between GB and MC, remarkably similar to the 
estimate of 0.01 that we found with approximately eight 
times more sampling effort but a subset of the microsatel- 



2013 The Authors. Ecology and Evolution published by John Wiley & Sons Ltd. 



3159 



Polar Bear Genetics in the Canadian Arctic 



L. Campagna ef al. 




-50° -30° -10° 10° 30° 



Figure 4. Phylogeographic analysis of polar bears using CR sequences. (A) Maximum parsimony network showing 95% probability linkages 
among 22 CR haplotypes obtained from 131 individuals sampled across the species' range, including sequences from Edwards et al. (2011). 
Localities were grouped according to geographic proximity (for simplicity in displaying the results) and color coded in the network. Other details 
as in Figure 3A. (B) Bayesian phylogenetic topology with posterior probabilities indicating node support. Individuals are color coded by locality as 
in Figure 4A. (C) Map showing sampling localities. 



Table 1. Demographic population parameters estimated using IMa2; population size parameters (6 = 4W/i, where N is the effective population 
size) for GB, MC, and the ancestral population, splitting time multiplied by the mutation rate (f/i), and effective number of migrants per genera- 
tion (2Wm). 







"mc 


^Ancestral 


tp 


2Nm 


mtDNA 


0.38 [0^1.43] 


2.08 [0.08-38.08] 


10.73 1 


0.14 1 


1.59 [0.12-29.62] 


nuDNA sexes 


2.98 [0.43-12.18] 


0.93 [0.13-4.18] 


7.28 [4.93-11.18] 


1.72 x 10~ 3 [2.8 x 10~ 4 -7 x 10~ 3 ] 


0.12 [0-19.7] 


combined 












nuDNA females 


1.88 [0.13-13.57] 


0.33 [0.08-2.18] 


6.48 [4.53-8.68] 


2.8 x 10~ 3 [1.2 x 10~ 4 -0.01] 


1.91 [0-15.69] 


nuDNA males 


1.73 1 


0.33 1 


7.23 [5.08-9.43] 


4 x 10~ 5 [0-0.03] 


0.84 [0-20.40] 



Values where the posterior probability peaks as well as 95% highest posterior density intervals ([95% HPD]: the shortest interval that contains 
95% of the posterior probability) are shown. Parameters were calculated independently using mtDNA and nuDNA, and for the latter we obtained 
estimates for both sexes separately and combined. 

1 Posterior density reaches lower values but not 0 near the upper or lower limit of the prior. 
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lites (nine out of 16). The differentiation between our two 
focal populations is moderate compared to that found by 
Paetkau et al. (1999) over a broader spatial scale: pairwise 
F sr range 0.002-0.108. Census population sizes estimated 
using mark-recapture data were 1523 ± 285 individuals 
for GB and 284 ± 59 bears for MC (Furnell and Schwe- 
insburg 1984; Taylor et al. (2006, 2009); Obbard et al. 
2010). Frankham (1995) found effective population sizes 
(N e ) to be often markedly smaller than census population 
sizes (N), particularly for bears where NJN was between 
0.28 and 0.69. In a study conducted on polar bears from 
the Beaufort Sea (Alaska), Cronin et al. (2009) found NJ 
N to be -0.182. Similarly, our estimations of effective 
population size using approximate Bayesian computation 
suggest a NJN ratio of 0.28 for both GB and MC (NJN: 
421/1523 vs. 80/284). Thus, with relatively small effective 
population sizes (especially for MC where polar bear 
numbers were reduced drastically by overhunting) it is 
likely that genetic drift has played a role in shaping differ- 
ences in allele frequencies. Although the results from our 
IMa2 analysis must be interpreted carefully, the data sug- 
gest that the genetic differences between these populations 
arose within the last hundreds of years. This time frame 
suggests that the land barrier that separates GB and MC 
(the Boothia Peninsula) is sufficient to at least partially 
restrict polar bear movement across these areas, and that 
the Bellot Strait does not constitute a corridor connecting 
GB and MC. We agree with Paetkau et al. (1999) that it is 
surprising that an animal with such a striking capacity to 
travel long distances would show genetic evidence of 
restricted movements across relatively small geographic 
areas. As Paetkau et al. (1999) hypothesized, this is possi- 
bly a consequence of sea ice and prey distribution. Addi- 
tionally, the fidelity females show to denning areas 
(Amstrup and Gardner 1994; Zeyl et al. 2010) could play a 
role in shaping population differentiation. 

In their broad-scale study of the origin of modern polar 
bears, Edwards et al. (2011) found few highly supported 
nodes (posterior probability >0.85) within a CR clade 
comprising ~50 ancient and modern polar bears. Despite 
sampling individuals throughout the species' range (see 
Fig. 4C), their study did not find an unequivocal phyloge- 
ographic pattern. Both Edwards et al. (2011) and Hailer 
et al. (2012) have suggested that multiple hybridization 
events occurred between brown bears and polar bears 
throughout the Pleistocene, and that mitochondria have 
introgressed from the former to the latter species. Thus, it 
is possible that such mtDNA introgression confounded our 
attempts to infer polar bear demographic history. Using 
CR sequences from GB and MC individuals we found 
marked genetic differentiation between our inferred focal 
populations ($ sx = 0.47, F ST = 0.29), yet no clear phyloge- 
ographic pattern when compared to the CR sequences 



obtained by Edwards et al. (2011). This raises the interest- 
ing question of how such frequency-based differences in 
mtDNA between GB and MC arose when population struc- 
ture seems to be absent at a larger scale. Edwards et al. 
(2011) interpret the distribution of modern bear matrilines 
as the product of long periods of stability in local condi- 
tions interrupted by rapid dispersal events and occasional 
hybridization with brown bears. Thus, the differences we 
observed in our mitochondrial data could have originated 
as populations became isolated in different glacial refugia, 
as has been suggested to explain mitochondrial genetic pat- 
terns in other Arctic species (e.g., rock ptarmigans, Holder 
et al. 1999; collared lemmings, Federov and Stenseth 2002). 
A possible explanation for the mitochondrial differentia- 
tion apparent between GB and MC populations is that it is 
the product of a founder effect. These areas could have 
been colonized by small groups of individuals possessing 
mitochondrial haplotypes that diverged elsewhere; fine-scale 
differentiation was produced by random sorting of mtDNA 
diversity that existed at the time of colonization. Phylogeo- 
graphic studies of other Arctic species suggest modern 
ranges are the product of postglacial expansions from mul- 
tiple glacial refugia where populations had remained iso- 
lated (e.g., Holder et al. 1999; Federov and Stenseth 2002). 
Ultimately, more intensive geographic sampling across as 
many polar bear populations as possible, combined with 
phylogeographic data from an array of species with diver- 
gent life histories (especially less vagile species) will help us 
understand the dynamic history of polar bear populations. 

Four lines of evidence suggest that the mitochondrial 
and nuclear genetic patterns observed between GB and MC 
differ in origin. First, our IMa2 estimations of effective 
population sizes inferred through mtDNA and nuDNA are 
contradictory. While mtDNA suggests the effective popula- 
tion size of MC is approximately fivefold larger than that 
of GB, nuDNA shows a threefold difference in the opposite 
direction (there is a fivefold difference in favor of GB when 
N e is estimated with ONeSAMP using nuDNA). Direct esti- 
mates of current census population sizes in the study area 
(Furnell and Schweinsburg 1984; Taylor et al. (2006, 2009); 
Obbard et al. 2010) indicate that in fact the population of 
GB is roughly fivefold larger than that of MC. Assuming 
effective population sizes are proportional to the absolute 
number of individuals in the population (although gener- 
ally lower; Frankham 1995), then this is consistent with the 
estimates obtained using nuDNA that show larger effective 
population sizes in GB, supporting our interpretation that 
mtDNA reveals differentiation at a different time scale 
within the species. Second, while our IMa2 analysis using 
mtDNA shows that the GB and MC populations have expe- 
rienced gene flow since their establishment, the nuclear 
data imply that the contrary is true. A third line of evidence 
comes from the geographic scale at which nuDNA and 
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mtDNA genetic divergence is observed. While DNA micro- 
satellites reveal differentiation at a local scale (e.g., the 
study by Paetkau et al. (1999) shows higher F sr values 
when comparing geographically more distant populations 
consistent with population-level isolation by distance sensu 
Wright 1943), the differentiation in mtDNA is common to 
individuals across the species' range. Finally, our IMa2 esti- 
mates of divergence times are also discordant. Although 
the limited differentiation of nuDNA suggests that these 
populations were established during the last century, the 
divergence in mtDNA probably dates to before or during 
the Holocene. Edwards et al. (2011) dated the matrilineal 
ancestor of all extant polar bears (not only those from MC 
and GB) to 51,000-20,000 years before present; this coales- 
cence time implies that the genetic diversity observed in 
our CR sequences may be much older than we found, per- 
haps consistent with the long tail in the splitting time dis- 
tribution. Again, our estimates of divergence times are 
approximate as IMa2 assumes that all populations that are 
exchanging genes have been sampled and this is clearly not 
true of our surveys. Our analyses did not include all polar 
bear populations interacting in the model at the same time; 
however, conducting analyses across the species' range with 
the sampling intensity we have achieved for MC and GB 
would be both logistically challenging and computationally 
intense. Strasburg and Rieseberg (2010) found IMa is 
robust to small or moderate violations of the IM model 
assumptions, including introgression from an unsampled 
population. Nevertheless, future effort should focus on val- 
idating our findings by expanding sampling to the neigh- 
boring populations found by Paetkau et al. (1999) to be 
most closely related to GB and MC, and thus likely to be 
exchanging genes with them. In summary, our findings 
illustrate how nuclear and mitochondrial DNA markers 
may capture different aspects of the demographic history 
of a species. Thus, results from analyses of either class of 
marker considered alone may only reveal part of a species 
demographic and evolutionary history. Moreover, the 
combination of these markers in the same analysis must be 
interpreted with caution as they could potentially produce 
intermediate spurious results. 

Gene flow between GB and MC populations 

Our indirect estimates of gene flow using IMa2 suggest 
gene exchange has been low since the GB and MC popula- 
tions began to differentiate. However, the effective number 
of migrants per generation (2Nm) calculated using CR 
sequences was significantly different from zero. As we 
argued above, we do not think that this represents the sig- 
nature of contemporary gene flow between GB and MC, 
but rather the vestige of genealogical events that occurred 
elsewhere. Our estimates of the effective number of 



migrants per generation using nuDNA did not differ signif- 
icantly from zero suggesting that contemporary gene flow 
between these populations is insufficient to homogenize 
them. This implies that low microsatellite differentiation 
and the admixture shown by the Bayesian assignment anal- 
yses (Fig. 2) is a product of recent separation rather than 
ongoing gene flow. Interestingly, there were three cases in 
which adult females initially caught in MC were recaptured 
in GB, two of which had a pair of cubs. We do not know if 
these represent real dispersal events or simply transitory 
movements of individuals between populations that did 
not lead to gene flow. Taylor et al. (2001) reported limited 
movement among populations north of our study area. For 
example, they found that only two of 65 bears for which 
there were multiyear satellite telemetry data available 
moved between defined populations among years. Annual 
rates of movement across population boundaries based on 
mark-recapture data were estimated at between 0.4% and 
8.9%. Certainly, if the movement between MC and GB that 
we observed led to gene flow, their genetic contribution to 
the overall population was not sufficient enough to be evi- 
dent in our IMa2 analysis. 

Finally, we found no compelling evidence for sex-biased 
dispersal or gene flow. Our mark-recapture data suggest 
that males and females do not differ in the distances trav- 
eled between years. The IMa2 results did not show evi- 
dence of gene flow between populations either for analyses 
with sexes combined or when analyzed separately. These 
results are consonant with those of Taylor et al. (2001) 
who did not find differences in distances moved between 
sexes in six more northerly Arctic polar bear populations. 
Zeyl et al. (2009) combined microsatellite and field data 
to investigate spatial patterns of relatedness and dispersal 
in polar bears of the Barents Sea. Their analysis of genetic 
relatedness of same-sex dyads revealed stronger kin struc- 
ture in female dyads compared with male dyads. Zeyl 
et al. (2009) concluded that at broader spatial scales effec- 
tive dispersal was slightly male biased, and thus that male- 
biased gene flow may preclude genetic differentiation. We 
found that females had stronger signal of nuDNA differen- 
tiation between MC and GB, both in our Bayesian cluster- 
ing analyses and reflected in females showing larger 
estimates of tfi in the IMa2 analyses. This result would be 
consistent with male-biased gene flow or could simply be 
spurious, reflecting the lower sample sizes for males. 

Conservation and management implications 

Both nuDNA and mtDNA show that polar bears from 
the Canadian Arctic Archipelago spanning GB and MC 
comprise two genetically distinct populations. The level 
of differentiation between GB and MC is moderate when 
compared to polar bears from other populations (Pae- 
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tkau et al. 1999). However, our data suggest that gene 
flow is currently insufficient to homogenize these popula- 
tions. These genetic differences imply that GB and MC 
should be treated as separate management units (sensu 
Moritz 1994) supporting previous findings (Taylor and 
Lee 1995; Bethke et al. 1996; Paetkau et al. 1999; Taylor 
et al. 2001). Our mark-recapture data suggest that indi- 
viduals from MC travel on average larger distances than 
those from GB, although statistical significance was only 
achieved when the three females that moved between 
populations were included. This could reflect differences 
in the availability of high quality habitat for seals (Barber 
and Iacozza 2004) or fewer mating opportunities, the 
later being a consequence of the larger area and lower 
polar bear density of MC compared to GB (Obbard et al. 
2010). Current efforts to monitor and model the demo- 
graphics of GB and MC separately are justified and cru- 
cial to the recovery of the MC population from intense 
hunting. Future effort should include the neighboring 
populations in a more comprehensive genetic and ecolog- 
ical analysis that would provide a better understanding of 
the connectivity of polar bear populations in the Cana- 
dian Arctic Archipelago. 

Acknowledgments 

For assistance in the field, we are indebted to Joe Ashe- 
vak, Julie Beauchesne, Jacob Keanik, Remi Krikort, Joani 
Kringayark (deceased), Peter Krizan, Jody M'Crae, Steve 
Millar, Brian Morrison, Susan Polischuk, Toni Romito, 
and Mitch Taylor. We thank Miranda Hunt and Candace 
Scott for providing technical expertise in the lab. This 
study was funded by the Department of Environment of 
the Government of Nunavut (formerly the Department of 
Sustainable Development), the Department of Environ- 
ment and Natural Resources of the Government of the 
Northwestern Territories (formerly the Department of 
Renewable Resources), the Nunavut Wildlife Management 
Board, Indian and Northern Affairs (Canada), the Univer- 
sity of Saskatchewan, and by NSERC Discovery Grants to 
P. T. B. and S. C. L. We thank the Polar Continental 
Shelf Program (Natural Resources Canada) for Logistical 
Support (in-kind funding). L. C. thanks the Senate Advi- 
sory Research Committee at Queen's University for post- 
doctoral funding and B. L. S. thanks Queen's University 
for Graduate Studies Awards. We thank Matthew A. 
Cronin and an anonymous reviewer for their helpful 
comments on a previous version of this manuscript. 

Data Accessibility 

Our microsatellite data set and all Structure, Geneland, 
Arlequin, and IMa2 input files are available on the Dryad 



Repository (http://dx.doi.org/10.5061/dryad.hlv63). Mito- 
chondrial sequences were deposited in GenBank (Acces- 
sion Numbers KF192517 - KF 192602). 

Author Contribution 

L. C, P. J. v. C. d. G., P. T. B., and S. C. L. designed 
research; L. C, B. L. S., S. N. A., D. S. W., and M. G. D. 
performed research; L. C, P. J. v. C. d. G., and S. C. L. 
analyzed the data and wrote the paper with help from the 
other authors. 

Conflict of Interest 

None declared. 
References 

Amstrup, S. C, and C. Gardner. 1994. Polar bear maternity 

denning in the Beaufort Sea. J. Wildl. Manage. 58:1-10. 
Amstrup, S. C, B. G. Marcot, and D. C. Douglas. 2008. A 

Bayesian network modeling approach to forecasting the 21st 

century worldwide status of polar bears. Pp. 213-268 in E. 

T. DeWeaver, C. M. Blitz and L. B. Tremblay, eds. Arctic 

sea ice decline: observations, projections, mechanisms, and 

implications. Geophys. Monogr. Ser. 180. 
Barber, D. G., and J. Iacozza. 2004. Historical analysis of sea 

ice conditions in M'Clintock Channel and the Gulf of 

Boothia, Nunavut: implications for ringed seal and polar 

bear habitat. Arctic 57:1-14. 
Bethke, R., M. Taylor, S. Amstrup, and F. Messier. 1996. 

Population delineation of polar bears using satellite collar 

data. Ecol. Appl. 6:311-317. 
Clement, M., D. Posada, and K. A. Crandall. 2000. TCS: a 

computer program to estimate gene genealogies. MoL Ecol. 

9:1657-1659. 

Crompton, A. E., M. E. Obbard, S. D. Petersen, and 

P. J. Wilson. 2008. Population genetic structure in polar 
bears (Ursus maritimus) from Hudson Bay, Canada: 
implications for future climate change. Biol. Conserv. 
141:2528-2539. 

Cronin, M. A., S. C. Amstrup, and K. T. Scribner. 2006. 
Microsatellite DNA and mitochondrial DNA variation in 
polar bears in the Beaufort and Chukchi seas, Alaska. Can. 
I. Zool. 84:655-660. 

Cronin, M. A., S. C. Amstrup, S. L. Talbot, G. K. Sage, and 
K. S. Amstrup. 2009. Genetic variation, relatedness, and 
effective population size of polar bears (Ursus maritimus) 
in the southern Beaufort Sea, Alaska. J. Hered. 100:681- 
690. 

Derocher, A. E., N. I. Lunn, and I. Stirling. 2004. Polar bears 
in a warming climate. Integr. Comp. Biol. 44:163-176. 

Durner, G. M., D. C. Douglas, R. M. Nielson, S. C. Amstrup, 
T. L. McDonald, I. Stirling, et al. 2009. Predicting 



2013 The Authors. Ecology and Evolution published by John Wiley & Sons Ltd. 



3163 



Polar Bear Genetics in the Canadian Arctic 



L. Campagna ef al. 



21st-century polar bear habitat distribution from global 

climate models. Ecol. Monogr. 79:25-58. 
Earl, D. A., and B. M. vonHoldt. 2012. Structure harvester: a 

website and program for visualizing STRUCTURE output 

and implementing the Evanno method. Conserv. Genet. 

Resour. 4:359-361. 
Edwards, C. J., M. A. Suchard, P. Lemey, J. J. Welch, 

I. Barnes, T. L. Fulton, et al. 2011. Ancient hybridization 

and an Irish origin for the modern polar bear matriline. 

Curr. Biol. 21:1251-1258. 
Evanno, G., S. Regnaut, and J. Goudet. 2005. Detecting the 

number of clusters of individuals using the software 

STRUCTURE: a simulation study. Mol. Ecol. 14:2611-2620. 
Excoffier, L., and H. E. L. Lischer. 2010. Arlequin suite ver 3.5: 

a new series of programs to perform population genetics 

analyses under Linux and Windows. Mol. Ecol. Resour. 

10:564-567. 

Federov, V. B., and N. C. Stenseth. 2002. Multiple glacial 

refugia in the North American Arctic: inference from 

phylogeography of the collared lemming (Dicrostonyx 

groenlandicus) . Proc. Biol. Sci. 269:2071-2077. 
Frankham, R. 1995. Effective population size/adult population 

size ratios in wildlife: a review. Genet. Res. 66:95-107. 
Furnell, D. J., and R. E. Schweinsburg. 1984. Population 

dynamics of central Canadian Arctic island polar bears. J. 

Wildl. Manage. 48:722-728. 
Gelman, A., and D. B. Rubin. 1992. Inference form iterative 

simulation using multiple sequences. Stat. Sci. 7:457-472. 
Goudet, J. 1995. FSTAT (version 1.2): a computer program to 

calculate F-statistics. J. Hered. 86:485-486. 
Guillot, G., F. Mortier, and A. Estoup. 2005. Geneland: a 

program for landscape genetics. Mol. Ecol. Notes 5:712- 

715. 

Hailer, F., V. E. Kutschera, B. M. Hallstrom, D. Klassert, 
S. R. Fain, J. A. Leonard, et al. 2012. Nuclear genomic 
sequences reveal that polar bears are an old and distinct 
bear lineage. Science 336:344-347. 

Hall, T. A. 1999. BIOEDIT: a user-friendly biological sequence 
alignment editor and analysis program for Windows 95/98/ 
NT. Nucleic Acids Symp. Ser. 41:95-98. 

Hasegawa, M., H. Kishino, and T. A. Yano. 1985. Dating of 
the human-ape splitting by a molecular clock of 
mitochondrial DNA. J. Mol. Evol. 22:160-174. 

Hey, J. 2010. Isolation with migration models for more than 
two populations. Mol. Biol. Evol. 27:905-920. 

Hey, J., and R. Nielsen. 2004. Multilocus methods for 

estimating population sizes, migration rates and divergence 
time, with applications to the divergence of Drosophila 
pseudoobscura and D. persimilis. Genetics 167:747-760. 

Holder, K., R. Montgomerie, and V. L. Friesen. 1999. A test of 
the glacial refugium hypothesis using patterns of 
mitochondrial and nuclear DNA sequence variation in rock 
ptarmigan (Lagopus mutus). Evolution 53:1936-1950. 



Huang, Q. Y., F. H. Xu, H. Shen, H. Y. Deng, Y. J. Liu, Y. Z. Liu, 

et al. 2002. Mutation patterns at dinucleotide microsatellite 

loci in humans. Am. J. Hum. Genet. 70:625-634. 
Huelsenbeck, J. P., and F. Ronquist. 2001. MRBAYES: Bayesian 

inference of phylogenetic trees. Bioinformatics 17:754-755. 
Kalinowski, S. T., M. L. Taper, and T. C. Marshall. 2007. 

Revising how the computer program CERVUS 

accommodates genotyping error increases success in 

paternity assignment. Mol. Ecol. 16:1099-1106. 
Kerr, R. A. 2002. A warmer Arctic means change for all. 

Science 297:1490-1492. 
Langella, O. 1999. Populations 1.2.30: a population genetic 

software. CNRS UPR9034. Available at http://www.pge. 

cnrs-gif.fr/bioinfo/populations/ (accessed Februaiy 13, 2011). 
Lee, J., and M. Taylor. 1994. Aspects of the polar bear harvest 

in the Northwest Territories, Canada. Int. Conf. Bear Res. 

Manage. 9:237-243. 
Mantel, N. 1967. The detection of disease clustering and a 

generalized regression approach. Cancer Res. 27:209-220. 
Moritz, C. 1994. Defining 'evolutionary significant units' for 

conservation. Trends Ecol. Evol. 9:373-375. 
Ohta, T., and M. Kimura. 1973. The model of mutation 

appropriate to calculate the number of electrophoretically 

detectable alleles in a genetic population. Genet. Res. 

22:201-204. 

Paetkau, D., and C. Strobeck. 1994. Microsatellite analysis of 
genetic variation in black bear populations. Mol. Ecol. 
3:489^95. 

Paetkau, D., W. Calvert, I. Stirling, and C. Strobeck. 1995. 

Microsatellite analysis of population structure in Canadian 

polar bears. Mol. Ecol. 4:347-354. 
Paetkau, D., G. F. Shields, and C. Strobeck. 1998. Gene flow 

between insular, coastal and interior populations of brown 

bears in Alaska. Mol. Ecol. 7:1283-1292. 
Paetkau, D., S. C. Amstrup, E. W. Born, W. Calvert, 

A. E. Derocher, G. W. Garner, et al. 1999. Genetic structure 

of the world's polar bear populations. Mol. Ecol. 8:1571- 

1584. 

Peakall, R., and P. E. Smouse. 2006. GENALEX 6: genetic 

analysis in Excel. Population genetic software for teaching 

and research. Mol. Ecol. Notes 6:288-295. 
Obbard, M. E., G. W. Thiemann, E. Peacock, and T. D. 

DeBruyn. 2010. Polar bears: proceedings of the 15th 

working meeting of the IUCN/SSC polar bear specialist 

group, Copenhagen, Denmark, 29 June-3 July 2009. IUCN, 

Gland, Switzerland. 
Posada, D. 2008. jModelTest: phylogenetic model averaging. 

Mol. Biol. Evol. 25:1253-1256. 
Pritchard, J. K., M. Stephens, and P. Donnelly. 2000. Inference 

of population structure using multilocus genotype data. 

Genetics 155:945-959. 
Rice, W. R. 1989. Analyzing tables of statistical tests. Evolution 

43:223-225. 



3164 



© 201 3 The Authors. Ecology and Evolution published by John Wiley & Sons Ltd. 



L. Campagna ef al. 



Polar Bear Genetics in the Canadian Arctic 



Ronquist, F., and J. P. Huelsenbeck. 2003. MrBayes 3: Bayesian 
phylogenetic inference under mixed models. Bioinformatics 
19:1572-1574. 

Rousset, F. 2000. Genetic differentiation between individuals. J. 
Evol. Biol. 13:58-62. 

Saarma, U., S. Y. W. Ho, O. G. Pybus, M. Kaljuste, I. M. 
Tumanov, I. Kojola, et al. 2007. Mitogenetic structure of 
brown bears (Ursus arctos L.) in north-eastern Europe and a 
new time-frame for the formation of brown bear lineages. 
Mol. Ecol. 16:401-413. 

Saunders, B. L. 2005. The mating system of polar bears in the 
central Canadian Arctic. Dissertation (M.Sc), Queen's 
University, Kingston, Ontario. 

Schliebe, S., 0. Wiig, A. Derocher, and N. Lunn; IUCN SSC 
Polar Bear Specialist Group. 2008. Ursus maritimus in 
IUCN 2012. IUCN Red List of Threatened Species. Version 
2012.2. Available at www.iucnredlist.org (accessed January 
3, 2013). 

Shriver, M. D., L. Jin, E. Boerwinkle, R. Deka, R. E. Ferrell, 
and R. Chakraborty. 1995. A novel measure of genetic 
distance for highly polymorphic tandem repeat loci. Mol. 
Biol. Evol. 12:914-920. 

Stirling, I., and A. E. Derocher. 1993. Possible impacts of 
climatic warming on polar bears. Arctic 46:240-245. 

Strasburg, J. L., and L. H. Rieseberg. 2010. How robust 
are "isolation with migration" analyses to violations of 
the IM model? A simulation study. Mol. Biol. Evol. 
27:297-310. 

Taberlet, P., J.- J. Camarra, S. Griffin, E. Uhres, O. Hanotte, 
L. P. Waits, et al. 1997. Noninvasive genetic tracking of the 
endangered Pyrenean brown bear population. Mol. Ecol. 
6:869-876. 

Tallmon, D. A., A. Koyuk, G. Luikart, and M. A. Beaumont. 

2008. ONeSAMP: a program to estimate effective 

population size using approximate Bayesian computation. 

Mol. Ecol. Resour. 8:299-301. 
Taylor, M. K., and L. J. Lee. 1995. Distribution and abundance 

of Canadian polar bear populations: a management 

perspective. Arctic 48:147-154. 
Taylor, M. K., S. Akeeagok, D. Andriashek, W. Barbour, 

E. W. Born, W. Calvert, et al. 2001. Delineating Canadian 



and Greenland polar bear ( Ursus maritimus) populations by 
cluster analysis of movements. Can. J. Zool. 79:690-709. 

Taylor, M. K., J. Laake, P. D. McLoughlin, H. D. Cluff, and F. 
Messier. 2006. Demographic parameters and harvest-explicit 
population viability analysis for polar bears in M'Clintock 
Channel, Nunavut, Canada. J. Wildl. Manage. 70:1667-1673. 

Taylor, M. K., J. Laake, P. D. McLoughlin, H. D. Cluff, and F. 
Messier. 2009. Demography and population viability of 
polar bears in the Gulf of Boothia, Nunavut. Mar. Mamm. 
Sci. 25:778-796. 

Wilgenbusch, J. C, D. L. Warren, and D. L. Swofford. 2004. 
AWTY: a system for graphical exploration of MCMC 
convergence in Bayesian phylogenetic inference. Available at 
http://ceb.csit.fsu.edu/awty (accessed November 8, 2012). 

Wright, S. 1943. Isolation by distance. Genetics 28:114-138. 

Zeyl, E., J. Aars, D. Ehrich, and 0. Wiig. 2009. Families in 
space: relatedness in the Barents Sea population of polar 
bears (Ursus maritimus). Mol. Ecol. 18:735-749. 

Zeyl, E., D. Ehrich, J. Aars, L. Bachmann, and 0. Wiig. 2010. 
Denning-area fidelity and mitochondrial DNA diversity of 
female polar bears (Ursus maritimus) in the Barents Sea. 
Can. J. Zool. 88:1139-1148. 

Supporting Information 

Additional Supporting Information may be found in the 
online version of this article: 

Appendix SI. Plots of the likelihood of the data for dif- 
ferent values of K in the Structure analyses together with 
a table indicating the proportion of membership of each 
predefined population to the inferred clusters (when 
K= 2). 

Appendix S2. Posterior density curves for the different 
parameters estimated using mitochondrial (CR sequences) 
and nuclear (DNA microsatellites) data in IMa2. 
Appendix S3. Calculations of the splitting time parameter 
in years using information from IMa2 runs and approxi- 
mate molecular clock calibrations for CR sequences and 
dinucleotide DNA microsatellites. 



2013 The Authors. Ecology and Evolution published by John Wiley & Sons Ltd. 



3165 



